Spatial Correlation Analysis Between Conflict and Food Prices in The Sahel Region between 2010-2020

By Syver Jahren Petersen

2022-08-11

1. Introduction

In this report I conduct a study of the spatial relationship between food prices and conflict events in the Sahel region.

In low-income communities food prices often have great importance as food costs make up a significant portion of household spending. Food prices are considered an important indicator when considering food security often employed as a proxy metric for assessing food security. Various violent conflicts can affect food markets drastically; impeding food distribution, destroying crops and disturbing or interrupting other food producing activities; increasing price volatility and rising food prices.

Additionally, food insecurity produces an environment vulnerable to armed conflict, exacerbating already tense situations. When resources are dwindling and perspectives for better futures are scarce, recruitment for paramilitary and terrorist groups are on the rise. This feedback loop between food prices and the occurrence of armed conflicts is an important dynamic to be explored.

Loading key packages

suppressPackageStartupMessages(library(tidyverse))
suppressPackageStartupMessages(library(lubridate))
suppressPackageStartupMessages(library(sf))
suppressPackageStartupMessages(library(zoo))
library(ggplot2)
library(ggspatial)
suppressPackageStartupMessages(library(spatstat))
library(stars)
# Loading required package: abind

2. Data Pre-Processing: Food Prices

This point based data is collected by The World Food Programme and scourced from Humanitarian Data Exhange (https://data.humdata.org/). The food prices are recorded on at a large number of markets in the study region on a montly basis.

Niger_price <- read.csv("C:/Users/Laptop/Documents/Munster/5. Data to Knowledge/Food Prices SAHEL/Data/Price_Data/with_lat_long/Niger_prices_with_lat_long.csv")

Mali_price <- read.csv("C:/Users/Laptop/Documents/Munster/5. Data to Knowledge/Food Prices SAHEL/Data/Price_Data/with_lat_long/Mali_prices_with_lat_long.csv")

Burkina_price <- read.csv("C:/Users/Laptop/Documents/Munster/5. Data to Knowledge/Food Prices SAHEL/Data/Price_Data/with_lat_long/BurkinaF_prices_with_lat_long.csv")

Price columns are converted from character to numeric values and data columns are set to a date format ymd (using lubridate), and data is filtered by dates 2010-01-01 and 2021-12-31

Niger_prices1 <- Niger_price %>%
  mutate(X = row_number()) %>%
  mutate(price = as.numeric(price)) %>%
  mutate(date = ymd(date)) %>% filter(date >= as.Date("2010-01-01") & date <= as.Date("2020-12-31"))

Mali_prices1 <- Mali_price %>%
  mutate(X = row_number()) %>%
  mutate(price = as.numeric(price)) %>%
  mutate(date = dmy(date)) %>% filter(date >= as.Date("2010-01-01") & date <= as.Date("2020-12-31")) 

Burkina_prices1 <- Burkina_price %>%
  mutate(X = row_number()) %>%
  mutate(price = as.numeric(price)) %>%
  mutate(date = ymd(date)) %>% filter(date >= as.Date("2010-01-01") & date <= as.Date("2020-12-31")) %>% arrange(desc(date))

The price data is filtered to include only retail prices of the commodities Maize, Millet, Sorghum. These three grains are selected because they have a high degree of data coverage from 2010 to 2020, and are major food commodity types in the region.

# variables for filtering food types:
food <- c("Maize - Retail", "Millet - Retail", "Sorghum - Retail")
food1 <- c("Maize (white) - Retail", "Millet - Retail", "Sorghum (white) - Retail")


Niger_prices4 <- Niger_prices3 %>%  filter(cmname %in% food)

Mali_prices4 <- Mali_prices3 %>% filter(cmname %in% food)

Burkina_prices4 <- Burkina_prices3  %>% filter(cmname %in% food1)

As some of the monthly price records seem to be exact copies of preceeding records I suspect they have been arbitrary inherited to make up for missing data. I therefore treat these as missing data and remove all exact duplicates of other records within the same quarter, market and commodity.

Also, to make sure I have enough records in each to calculate later calculate a standard deviation, all records that are within quarters containing less than 3 records are removed, e.g. if there are only 2 records in for a single quarter, single market, and single commodity.

remove_nonunique_records <- function(data_frame) {
  data_frame2 <- data_frame %>%
  group_by(quarter, mktname, cmname)
data_frame3 <- data_frame2 %>% mutate(distinct_prices = n_distinct(price))
data_frame4 <- data_frame3[data_frame3$distinct_prices > 2, ] %>% ungroup()
{
  return(data_frame4)
}
}

Niger_prices5 <- remove_nonunique_records(Niger_prices4)
Mali_prices5 <- remove_nonunique_records(Mali_prices4)
Burkina_prices5 <- remove_nonunique_records(Burkina_prices4)

3. Design of Metrics for Food Prices:

In order to analyse the food prices meaningfully in relation to conflicts I had to aggregate them into a meaningful metric that represents how each market’s prices relates to the other markets prices. I chose two different metrics:

1) The mean standard deviation for all prices at each market within each quarter.

2) How each market per quarter deviates in percent from the mean of all prices within that quarter in that country.

First Food Price Metric

A standard deviation of prices are calculated quarterly for each commodity, at each market based on min-max normalized prices for each commodity.

calculate_standard_dev <- function(data_frame) {
  data_frame2 <- data_frame  %>%
  group_by(cmname) %>%
  mutate(price_norm = (price - min(price)) / (max(price) - min(price))) %>%
  ungroup()
  
  data_frame3 <- data_frame2 %>%
  group_by(quarter, mktname, cmname) %>% 
  mutate(stdev = sd(price_norm, na.rm =TRUE)) %>%
  ungroup()

  {
    return(data_frame3)
  }
}

Niger_prices6 <- calculate_standard_dev(Niger_prices5)
Mali_prices6 <- calculate_standard_dev(Mali_prices5)
Burkina_prices6 <- calculate_standard_dev(Burkina_prices5)

To aggregate all the different quarterly commodity standard deviations at each market I calculate the mean of the standard deviations for all the crops at each market, per quarter.

calculate_mean_stdev <- function(data_frame) {
  data_frame2 <- data_frame %>%
  group_by(quarter, mktname) %>%
  mutate(stdev_mean = mean(stdev, na.rm = TRUE)) %>% 
  ungroup()
  {
    return(data_frame2)
  }
}

Niger_prices7 <- calculate_mean_stdev(Niger_prices6)
Mali_prices7 <- calculate_mean_stdev(Mali_prices6)
Burkina_prices7 <- calculate_mean_stdev(Burkina_prices6)

To avoid duplicate locations when later doing a spatial join, all markets situated at the same location are merged.

merge_same_location_markets <- function(data_frame) {
data_frame2 <- data_frame %>%
  group_by(quarter, lat, long) %>% 
  mutate(stdev_mean = mean(stdev_mean, na.rm =TRUE)) %>%
  ungroup()

{
  return(data_frame2)
}
}

Niger_price8 <- Niger_prices7 #does not have any location markets

Mali_price8 <- merge_same_location_markets(Mali_prices7)

Burkina_price8 <- Burkina_prices7 #does not have any location markets

Second Food Price Metric

for each single quarter for each single market the percent difference from the mean of food prices per country per quarter is calculated.

mean_diff_column <- function(data_frame) { 
data_frame2<- data_frame %>% 
  group_by(quarter, cmname) %>% 
  mutate(avg_nat_mean= mean(price)) %>%
  ungroup()

data_frame3<- data_frame2%>% 
  group_by(quarter, cmname, mktname) %>% 
  mutate(avg_local_price= mean(price)) %>%
  ungroup()

data_frame4 <-data_frame3 %>% mutate(Percent_Difference = (avg_local_price - avg_nat_mean)/avg_nat_mean*100)

data_frame5<-data_frame4%>%
  group_by(quarter, mktname) %>%
  mutate(avg_difference = mean(Percent_Difference)) %>%
  ungroup()

{
  return(data_frame5)
}
}

Niger_price10 <- mean_diff_column(Niger_price8)
Mali_price10 <- mean_diff_column(Mali_price8)
Burkina_price10 <- mean_diff_column(Burkina_price8)

4. Data Pre-Processing: Conflict Events

This point based data is was collected The Armed Conflict Location & Event Data Project (ACLED) and sourced from their website (https://acleddata.com/). The conflict events are recorded based on reports in legacy media and social media content. Each record represents an event, however, if a battle occurred from January 10th to January 14th then the “battle” would have 1 unique record for each day.

All_conflicts <- read.csv("C:/Users/Laptop/Documents/Munster/5. Data to Knowledge/Food Prices SAHEL/Data/Conflict_Data/All_Conflict_Events.csv")

Niger_conflicts <- All_conflicts %>%
  select(c(data_id, event_date, event_type, country, latitude, longitude)) %>%
  filter(country == "Niger")

Burkina_conflicts <- All_conflicts %>%
  select(c(data_id, event_date, event_type, country, latitude, longitude)) %>%
  filter(country == "Burkina Faso")

Mali_conflicts <- All_conflicts %>%
  select(c(data_id, event_date, event_type, country, latitude, longitude)) %>%
  filter(country == "Mali")

To prepare the conflict data for further analysis I need to do some pre-processing such as create a quarter column in the data frame and arrange the data frame by date.

I also filter the conflict events to only include “violent conflicts”, eg. “Battles”, “Violence against civilians”, “Explosions/Remote violence”.


# add a "quarter" column
Niger_conflicts$quarter = Niger_conflicts$event_date
Mali_conflicts$quarter = Mali_conflicts$event_date
Burkina_conflicts$quarter = Burkina_conflicts$event_date

# function for pre-processing the conflict data
conflicts_filter = c("Battles", "Violence against civilians", "Explosions/Remote violence")
#conflicts_filter = c("Battles", "Protests", "Violence against civilians", "Explosions/Remote violence", "Riots")

Pre_process_conflicts <- function(data_frame) {
  
  data_frame2 <- data_frame %>% mutate(event_date = dmy(event_date)) %>%
  rename("date" = "event_date", "type" = "event_type", "lat" = "latitude", "long" = "longitude") %>% filter(type %in% conflicts_filter)
  
  data_frame2$date <- as.yearmon(data_frame2$date, "%Y%m")
  data_frame2$quarter <- as.yearqtr(dmy(data_frame2$quarter), "%Y%m")
  
  final_pre_processsed <- data_frame2 %>% arrange(date)
  
  {
    return(final_pre_processsed)
  }
}

# pre-process the conflict data
Niger_conflicts3 <- Pre_process_conflicts(Niger_conflicts)

Mali_conflicts3 <- Pre_process_conflicts(Mali_conflicts)

Burkina_conflicts3 <- Pre_process_conflicts(Burkina_conflicts)

5. Make Simple Feature Geometry Collections from Coordinates

The latitude and longitude columns of the conflict data and the food price data needs to be formatted into Well Known Text (WKT) format, which enables conversion into simple features geometry collections.

WKT_conversion_function <- function(data_frame) {
  data_frame$geometry <- paste("POINT(",data_frame$long, ",", data_frame$lat,")")

  data_frame2 <- data_frame %>% mutate(geometry = str_replace_all(geometry, " ", ""))

  data_frame3 <- data_frame2 %>% mutate(geometry = str_replace_all(geometry, ",", " "))  

  {
    return(data_frame3)
  }
}

WKTcol_to_sfc_and_make_sf_collection <- function(data_frame) {
  
  data_frame$geometry <- st_as_sfc(data_frame$geometry)
  final_sf_collection <- st_sf(data_frame, crs = "EPSG:4326")
  {
    return(final_sf_collection)
  }
}

Niger_conflicts4 <- WKT_conversion_function(Niger_conflicts3)
Niger_conflicts5 <- WKTcol_to_sfc_and_make_sf_collection(Niger_conflicts4)

Mali_conflicts4 <- WKT_conversion_function(Mali_conflicts3)
Mali_conflicts5 <- WKTcol_to_sfc_and_make_sf_collection(Mali_conflicts4)

Burkina_conflicts4 <- WKT_conversion_function(Burkina_conflicts3)
Burkina_conflicts5 <- WKTcol_to_sfc_and_make_sf_collection(Burkina_conflicts4)

Niger_price11 <- WKT_conversion_function(Niger_price10)
Niger_price12 <- WKTcol_to_sfc_and_make_sf_collection(Niger_price11)

Mali_price11 <- WKT_conversion_function(Mali_price10)
Mali_price12 <- WKTcol_to_sfc_and_make_sf_collection(Mali_price11)

Burkina_price11 <- WKT_conversion_function(Burkina_price10)
Burkina_price12 <- WKTcol_to_sfc_and_make_sf_collection(Burkina_price11)

6. Load Country Boundary Data

Niger_country <- st_read("C:/Users/Laptop/Documents/Munster/5. Data to Knowledge/Food Prices SAHEL/Data/Country_Boundaries/Niger_boundaries.geojsonl.json")
# Reading layer `Niger_boundaries.geojsonl' from data source 
#   `C:\Users\Laptop\Documents\Munster\5. Data to Knowledge\Food Prices SAHEL\Data\Country_Boundaries\Niger_boundaries.geojsonl.json' 
#   using driver `GeoJSON'
# Simple feature collection with 1 feature and 2 fields
# Geometry type: MULTIPOLYGON
# Dimension:     XY
# Bounding box:  xmin: 0.1529411 ymin: 11.69577 xmax: 15.97032 ymax: 23.51735
# Geodetic CRS:  WGS 84
Mali_country <- st_read("/Users/Laptop/Documents/Munster/5. Data to Knowledge/Food Prices SAHEL/Data/Country_Boundaries/Mali_boundaries.geojsonl.json")
# Reading layer `Mali_boundaries.geojsonl' from data source 
#   `C:\Users\Laptop\Documents\Munster\5. Data to Knowledge\Food Prices SAHEL\Data\Country_Boundaries\Mali_boundaries.geojsonl.json' 
#   using driver `GeoJSON'
# Simple feature collection with 1 feature and 2 fields
# Geometry type: MULTIPOLYGON
# Dimension:     XY
# Bounding box:  xmin: -12.26413 ymin: 10.14005 xmax: 4.235638 ymax: 24.99506
# Geodetic CRS:  WGS 84
Burkina_country <- st_read("C:/Users/Laptop/Documents/Munster/5. Data to Knowledge/Food Prices SAHEL/Data/Country_Boundaries/Burkina_boundaries.geojsonl.json")
# Reading layer `Burkina_boundaries.geojsonl' from data source 
#   `C:\Users\Laptop\Documents\Munster\5. Data to Knowledge\Food Prices SAHEL\Data\Country_Boundaries\Burkina_boundaries.geojsonl.json' 
#   using driver `GeoJSON'
# Simple feature collection with 1 feature and 2 fields
# Geometry type: MULTIPOLYGON
# Dimension:     XY
# Bounding box:  xmin: -5.522578 ymin: 9.391883 xmax: 2.390169 ymax: 15.07991
# Geodetic CRS:  WGS 84
All_countries <- st_read("C:/Users/Laptop/Documents/Munster/5. Data to Knowledge/Food Prices SAHEL/Data/Country_Boundaries/All_countries_boundaries.geojsonl.json")
# Reading layer `All_countries_boundaries.geojsonl' from data source 
#   `C:\Users\Laptop\Documents\Munster\5. Data to Knowledge\Food Prices SAHEL\Data\Country_Boundaries\All_countries_boundaries.geojsonl.json' 
#   using driver `GeoJSONSeq'
# Simple feature collection with 3 features and 2 fields
# Geometry type: MULTIPOLYGON
# Dimension:     XY
# Bounding box:  xmin: -12.26413 ymin: 9.391883 xmax: 15.97032 ymax: 24.99506
# Geodetic CRS:  WGS 84

7. Exploratory Data Analysis

Let’s explore the data I have so far. First, let’s look at the countries of interest:

  
ggplot() + 
  
  geom_sf(data = All_countries, fill= "antiquewhite") +
  annotation_scale(location = "br", width_hint = 0.5) + 
  annotation_north_arrow(location = "tl", which_north = "true", pad_x = unit(0.75, "in"), pad_y = unit(0.5, "in"), style = north_arrow_fancy_orienteering) +
  xlab("Longitude") + 
  ylab("Latitude") + 
  ggtitle("Countries of Interest") + 
  theme(plot.title = element_text(size=20, face="bold")) +
  theme(panel.grid.major = element_line(color = gray(.5), linetype = "dashed", size = 0.5), panel.background = element_rect(fill = "aliceblue")) +
  
  annotate(geom = "text", x = -2, y = 18, label = "Mali", fontface = "bold", color = "grey22", size = 4) +
  annotate(geom = "text", x = 10, y = 17, label = "Niger", fontface = "bold",color = "grey22", size = 4) +
  annotate(geom = "text", x = -1, y = 12, label = "Burkina Faso", fontface = "bold", color = "grey22", size = 4)
# Scale on map varies by more than 10%, scale bar may be inaccurate

The area I want to explore is the Sahel region, however, as this contains a wide array of countries with a high diversity of conflicts, economies, politics and culture, to make the study feasible, I chose to focus on only three countries in the western part of the Sahel, Mali, Niger and Burkina Faso.

As you can see the three countries are contiguous and share long borders with each other which is beneficial in terms of being able to treat it as a single geographic region and include cross-border effects between the two variables.

Now, let’s have a look at all the spatial distribution of violent conflict events between 2010-2020:

ggplot()  +
  geom_sf(data = All_countries, fill= "antiquewhite") +
  geom_sf(data = rbind(Niger_conflicts5, Mali_conflicts5, Burkina_conflicts5), size=0.3, aes(col = "red")) +
  annotation_scale(location = "br", width_hint = 0.5) + 
  annotation_north_arrow(location = "tl", which_north = "true", pad_x = unit(0.75, "in"), pad_y = unit(0.5, "in"), style = north_arrow_fancy_orienteering) +
  xlab("Longitude") + 
  ylab("Latitude") + 
  ggtitle("Spatial Distribution of Conflicts (2010-2020)") + 
  theme(plot.title = element_text(size=20, face="bold")) +
  theme(panel.grid.major = element_line(color = gray(.5), linetype = "dashed", size = 0.5), panel.background = element_rect(fill = "aliceblue")) +
  
  scale_color_manual(values=c("red"), labels=c("Conflicts")) + labs(color='Legend') 
# Scale on map varies by more than 10%, scale bar may be inaccurate

And the spatial distribution of the markets where the food prices have been recorded between 2010-2020:

ggplot() + 
  geom_sf(data = All_countries, fill= "antiquewhite") +
  geom_sf(data = rbind(Niger_price12, Mali_price12, Burkina_price12), size=1, aes(col = "navyblue"))+
  annotation_scale(location = "br", width_hint = 0.5) + 
  annotation_north_arrow(location = "tl", which_north = "true", pad_x = unit(0.75, "in"), pad_y = unit(0.5, "in"), style = north_arrow_fancy_orienteering) +
  xlab("Longitude") + 
  ylab("Latitude") + 
  ggtitle("Spatial Distribution of Food Markets (2010-2020)") + 
  theme(plot.title = element_text(size=20, face="bold")) +
  theme(panel.grid.major = element_line(color = gray(.5), linetype = "dashed", size = 0.5), panel.background = element_rect(fill = "aliceblue")) +
  scale_color_manual(values=c("navyblue"), labels=c("Markets")) + labs(color='Legend')
# Scale on map varies by more than 10%, scale bar may be inaccurate

Let’s have a look at the annual standard deviation changes:

Niger_mean_quarter <- Niger_price12 %>% group_by(quarter) %>% summarise(stdev_mean = mean(stdev_mean))
Mali_mean_quarter <- Mali_price12 %>% group_by(quarter) %>% summarise(stdev_mean = mean(stdev_mean))
Burkina_mean_quarter <- Burkina_price12 %>% group_by(quarter) %>% summarise(stdev_mean = mean(stdev_mean))
library("plotly")
# 
# Attaching package: 'plotly'
# The following object is masked from 'package:ggplot2':
# 
#     last_plot
# The following object is masked from 'package:stats':
# 
#     filter
# The following object is masked from 'package:graphics':
# 
#     layout

ggplotly(
ggplot() + 
  geom_line(data= Niger_mean_quarter, aes(y=stdev_mean, x=quarter, col="Niger")) + 
  
  geom_line(data=Mali_mean_quarter, aes(y = stdev_mean, x=quarter,col = "Mali")) + 
  
  geom_line(data= Burkina_mean_quarter, aes(y = stdev_mean, x=quarter,col="Burkina Faso")) + 

  labs(title="Quarterly Price Standard Deviation per Country", y="standard deviation of min-max normalised prices", x="quarter", 
       color=NULL) + scale_color_manual(labels = c("Niger", "Burkina Faso", "Mali"), 
                     values = c("Niger"="#00ba38", "Burkina Faso"="#f8766d", "Mali" = "blue")) +
  theme(axis.text.x = element_text( vjust=0.5 , size = 8), plot.title = element_text(size=16, face="bold"))
)

And, the quarterly conflicts per country:

Niger_quarterly_conflicts2 <- Niger_conflicts3 %>% group_by(quarter) %>% summarise(n_conflicts = n_distinct(data_id))
Mali_quarterly_conflicts2 <- Mali_conflicts3 %>% group_by(quarter) %>% summarise(n_conflicts = n_distinct(data_id))
Burkina_quarterly_conflicts2 <- Burkina_conflicts3 %>% group_by(quarter) %>% summarise(n_conflicts = n_distinct(data_id))


ggplotly(
ggplot() + 
  geom_line(data = Niger_quarterly_conflicts2, aes(x = quarter, y = n_conflicts, col="Niger")) +
  geom_line(data = Burkina_quarterly_conflicts2, aes(x = quarter, y = n_conflicts, col="Burkina Faso")) +
  geom_line(data = Mali_quarterly_conflicts2, aes(x = quarter, y = n_conflicts, col="Mali")) +
  
  labs(title="Quarterly Conflicts per Country", y="number of conflicts", x="quarter", color=NULL) + 
  scale_color_manual(labels = c("Niger", "Burkina Faso", "Mali"), 
                     values = c("Niger"="#00ba38", "Burkina Faso"="#f8766d", "Mali" = "blue")) +
  theme(axis.text.x = element_text( vjust=0.5 , size = 8), plot.title = element_text(size=16, face="bold"))
)

The annual overview shows that conflicts have risen somehow in the last 3-4 years. While the variation in food prices shows a less clear pattern with higher average price deviation between 2015-2020 than 2010-2015.

9. Analysis

Step 1: Kernel Point Density Surfaces

In order to relate the location of the conflict events to the food markets I calculated, based on the point data of the conflict events, kernel density surfaces for each quarter using Spatstat. A kernel density surface is a surface of values denoting the per area unit density of points over a restricted surface. In this case the outer geographic delineation (or “bounding box” if you want) is the outer border of the three countries merged together into one polygon. Thus, the density surfaces covered all the three countries which allowed for conflict events to have influence across countries’ international borders. The key parameter when calculating a kernel point density surface is the bandwidth of the kernel. Normally you would use an cross-validation based algorithm to optimize the bandwith for a given surface. However, as I wanted to calculate one surface for each quarter where each quarter had a different location and amount of points, each surface would have a different optimal bandwidth. I therefore ran the density calculations at different bandwidths between 5 km and 200 km. I use the default gaussian kernel for the density estimations.

Step 2: Extract Values from Point Density Surface at Market Locations

After calculating the point density surface from the conflict events, I extracted the density values at each market location and thus geographically linked the food prices and the conflict events together.

Step 3: Regression analysis

In the final step I ran a regression analysis between the conflict density value extracted from the point density surface, and the two food price metrics. The final output of the analysis was a regression coefficient between 1) the conflict event density values and the first price metric, and 2) the conflict event density values and the second price metric.

FINAL_results <- data.frame(mktname = 0, final_std = 0, v = 0, final_avg_diff = 0, geometry = 0)

merged_list_margets <- rbind(Niger_price12, Mali_price12, Burkina_price12)

merged_list_margets2 <- merged_list_margets %>% filter(quarter != "2021 Q1")

unique_quarters <- c(unique(merged_list_margets2$quarter))

# All_countries polygons
All_countries2 <- All_countries %>% st_transform(3857) %>% select(ADMIN)

# Bounding Box for All_countries3
bb <- st_geometry(All_countries2) %>% st_union()

BW_df <- data.frame()


# list of bandwidths for the point density maps
sigma_list <- list(200000, 150000, 100000, 90000, 80000, 70000, 60000 ,50000, 45000, 40000 ,35000, 30000, 25000, 20000, 15000, 10000, 5000)
# 
# dataframes for final correlation results per bandwidth  

sd_correlations <- data.frame(bandwidths = c(200000, 150000, 100000, 90000, 80000, 70000, 60000 ,50000, 45000, 40000 ,35000, 30000, 25000, 20000, 15000, 10000,5000), correlation = NA)

diff_from_mean_correlations <- data.frame(bandwidths = c(200000, 150000, 100000, 90000, 80000, 70000, 60000 ,50000, 45000, 40000 ,35000, 30000, 25000, 20000, 15000, 10000, 5000), correlation = NA)


#This loop runs all the steps of the analysis one time for each of the bandwidths at a time. For each bandwidth the analysis is also run one time for each quarter between 2010-2020
for (bandwidht in sigma_list) {

  for (Q in unique_quarters) {

# Step 1: Kernel Point Density Surfaces
    
  # Point data for conflicts for each country
  Niger_points <- Niger_conflicts5 %>% st_transform(3857) %>% filter(quarter == Q)
  Mali_points <- Mali_conflicts5 %>% st_transform(3857) %>% filter(quarter == Q)
  Burkina_points <- Burkina_conflicts5 %>% st_transform(3857) %>% filter(quarter == Q)

  # Format as point pattern data
  point_patern = c(bb, st_geometry(Niger_points), st_geometry(Mali_points), st_geometry(Burkina_points)) %>% as.ppp()


  skip_to_next <- FALSE
  
 # The density function takes the point pattern data and the bandwidth as arguments. The dimyx refers the pixel resolution
    tryCatch(density_map <- density(point_patern, sigma = bandwidht, dimyx = 1000, verbose = FALSE), error = function(e) { skip_to_next <<- TRUE})

    if(skip_to_next) { next }

  density_map_stars <- st_as_stars(density_map) %>% st_set_crs(3857)

  Niger_price13 <- Niger_price12 %>% filter(quarter == Q)
  Mali_price13 <- Mali_price12 %>% filter(quarter == Q)
  Burkina_price13 <- Burkina_price12 %>% filter(quarter == Q)

  # Market points
  All_markets_final <- rbind(Niger_price13, Mali_price13, Burkina_price13) %>% st_transform(3857)
  
  # count column
  All_markets_final$counts <- 1

  # Aggregate std of all commodities to one row per mktname/quarter:
  All_markets_final2 <- All_markets_final %>% group_by(mktname) %>% summarise(final_std = sum(stdev_mean/sum(counts)), final_avg_diff = sum(avg_difference/sum(counts)))

# Step 2: Extract values for each market point from density map:
  extracted <- st_extract(density_map_stars, All_markets_final2)

  Final_join <- st_join(All_markets_final2, extracted, join = st_equals)


  FINAL_results <- rbind(FINAL_results, Final_join)
  }
  
# Step 3: Regression analysis. For each bandwidth a correlation coefficient is calculated between the two metrics and teh conflict density extracted at the markets within each of the quarters. 
  sd_correlations[sd_correlations$bandwidths == bandwidht, "correlation"] <- with(na.omit(FINAL_results), cor(final_std, v))
  diff_from_mean_correlations[sd_correlations$bandwidths == bandwidht, "correlation"] <- with(na.omit(FINAL_results), cor(final_avg_diff, v))
  
  # The data frame is reset before running next bandwidth
  FINAL_results <- data.frame(mktname = 0, final_std = 0, v = 0, final_avg_diff = 0, geometry = 0)
}

PS: Some of the conflict points fall outside the country boundaries used as a bounding box and are therefore not included in the point density estimations. Warnings about duplicated points are ignored as this just means some conflict events have been recorded at the same location. All the points are still included in the density estimation.

As you can see in these three examples, each bandwidth gives a different density estimation from which the points are extracted (green points). As menstioned this is done for each quarter between 2010-2020.

Niger_price16 <- Niger_price12 %>% filter(quarter == "2016 Q1")
Mali_price16 <- Mali_price12 %>% filter(quarter == "2016 Q1")
Burkina_price16 <- Burkina_price12 %>% filter(quarter == "2016 Q1")

All_markets_final <- rbind(Niger_price16, Mali_price16, Burkina_price16) %>% st_transform(3857)

plot_bandwidths <- list(200000, 45000, 5000) 
# 

for (plot_bandwidth in plot_bandwidths) {

# Point data for conflicts for each country
  Niger_points <- Niger_conflicts5 %>% st_transform(3857) %>% filter(quarter == "2016 Q1")
  Mali_points <- Mali_conflicts5 %>% st_transform(3857) %>% filter(quarter == "2016 Q1")
  Burkina_points <- Burkina_conflicts5 %>% st_transform(3857) %>% filter(quarter == "2016 Q1")

  # Format as point pattern data
  point_patern = c(bb, st_geometry(Niger_points), st_geometry(Mali_points), st_geometry(Burkina_points)) %>% as.ppp()
  
  plot_density_map <- density(point_patern, sigma = plot_bandwidth, dimyx = 1000)
  
  plot(plot_density_map, main= paste("Point density map for Q1 2017 with bandwidth", plot_bandwidth), reset = FALSE)

  plot(All_markets_final, size=1, col = "green", add = TRUE)

}

9. Final Results

To visualize the tendency of the relationships between the final correlation coefficients and each used bandwidth I created two scatter plots with a smoothing lines with confidence bands.

ggplotly(
ggplot(sd_correlations, aes(x=bandwidths, y= correlation)) +
        geom_smooth()+
        geom_point() +
        ggtitle("Price Standard Deviation vs. Conflict Events")+
        ylab("Correlation coefficients") +
        xlab("Bandwidhts (meters)") +
        theme(plot.title = element_text(size=16, face="bold"))
)
ggplotly( ggplot(diff_from_mean_correlations, aes(x=bandwidths, y= correlation)) + geom_smooth() + geom_point() + ggtitle("Price Diff from National Mean vs. Conflict Events")+ ylab("Correlation coefficients") + xlab("Bandwidhts (meters)") + theme(plot.title = element_text(size=16, face="bold")) )

10. Discussion

The results show at smaller bandwidth distances, between 5 and 35 kilometers there seems to be a correlation between the market price difference from the national mean metric and the conflict event densities. As smaller bandwidths naturally include fewer points in each kernel, this would allude to that a high density of events within one quarter, happening close to a market, could cause prices to rise above the quarterly national mean.

On the other hand the correlation between the price variation or quarterly standard deviation of prices at each market and the conflict event densities shows an opposite relationship, higher correlation at larger bandwidth distances. This being said, the correlation coefficients are much lower than the ones of the deviation from mean metric, which makes this tendency less interesting and as it is more likely to be random.

The biggest uncertainty in the results comes from the fact that food prices are most of the time not recorded in markets in very close temporal and spatial proximity to significant violent conflict events.

Another limitation of the conflict event data is that it is based on reports from legacy media and social media content which means it is likely that many conflict events are not recorded as the flow of information is often limited from areas with ongoing violent conflict.

In conclusion, the market price difference from the national mean metric shows some promise in having a certain degree of explanatory power over the conflicts, or vice-versa. There is however, as discussed, great uncertainty in regards to how well the data is suited to conduct this analysis. To make stronger claims about these relationships more exploration is needed.